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In this paper, we propose an improvement of the adaptive biasing force (ABF) method, 
by projecting the estimated mean force onto a gradient. The associated stochastic process 
satisfies a non linear stochastic differential equation. Using entropy techniques, we prove 
exponential convergence to the stationary state of this stochastic process. We finally 
show on some numerical examples that the variance of the approximated mean force is 
reduced using this technique, which makes the algorithm more efficient than the standard 
ABF method. 
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1. Introduction 
1.1. The model 

Let us consider the Boltzmann-Gibbs measure : 

^i(dx) = Z~ 1 e~ l3v( ' x ^dx, (1.1) 

where x £ V N denotes the position of N particles in V. The space V is called the 
configuration space. One should think of I? as a subset of K™, or the n-dimensional 
torus T” (where T = K/Z denotes the one dimensional torus). The potential energy 
function V : V — > R. associates to the positions of the particles x € V its energy 

V{x). In addition, Z, x = f e~ pv ^ x) dx (assumed to be finite) is the normalization 

Jv 

constant and f3 = 1 /(ksT) is proportional to the inverse of the temperature T, ks 
being the Boltzmann constant. 


1 
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The probability measure p, is the equilibrium measure sampled by the particles 
in the canonical statistical ensemble. A typical dynamics that can be used to sample 
this measure is the Overdamped Langevin Dynamics: 

dX t = -W(A : t )dt + dWt , (1.2) 

where X t £ V N and Wt is a An-dimensional standard Brownian motion. Under 
loose assumptions on V, the dynamics ( X t )t>o is ergodic with respect to the measure 
/i, which means: for any smooth test function ip, 

lim l- [ ip(X t )dt = [ ipdfi, (1.3) 

T^+oo 1 J Q J 

i.e. trajectory averages converge to canonical averages. 


1.2. Metastability, reaction coordinate and free energy 


In many cases of interest, there exists regions of the configuration space where 


the dynamics (1.2) remains trapped for a long time, and jumps only occasionally 


to another region, where it again remains trapped for a long time. This typically 
occurs when there exist high probability regions separated by very low probability 
areas. The regions where the process (X t ) t >o remains trapped for very long times, 
are called metastable. 


Because of the metastability, trajectorial averages (1.3) converge very slowly to 


their ergodic limit. Many methods have been proposed to overcome this difficulty, 
and we concentrate here on the Adaptive Biasing Force (denoted ABF) method 
(see mmy In order to introduce the ABF method, we need another ingredient: a 
reaction coordinate (also known as an order parameter), £ = (£ 1; ..., £ m ) : D —> M m , 
£( 2 ) = z, where m < nN. Typically, in (1.2), the time-scale for the dynamics on 


£(A' t ) is larger than the time-scale for the dynamics on X t due to the metastable 
states, so that £ can be understood as a function such that £(A t ) is in some sense 
a slow variable compared to X t . We can say that £ describes the metastable states 
of the dynamics associated to the potential V. For a given configuration x, £(x) 
represents some macroscopic information. For example, it could represent angles or 
bond lengths in a protein, positions of defects in a material, etc ... In any case, it is 
meant to be a function with values in a small dimensional space (i.e. m < 4), since 
otherwise, it is difficult to approximate accurately the associated free energy which 


is a scalar function defined on the range of £ (see equation (1.5) below). The choice 
of a ’’good” reaction coordinate is a highly debatable subject in the literature. One 
aim of the mathematical analysis conducted here or in previous papers (see for 
example i) is to quantify the efficiency of numerical algorithms once a reaction 
coordinate has been chosen. 

The image of the measure p by £ is defined by: 


£*//:= exp(— (3A(z))dz, 


(1.4) 












January 29, 2015 1:49 PABF2D 


Convergence of an ABF method: Variance reduction by Helmholtz projection 3 


where A is the so-called free energy associated with the reaction coordinate £. By 
the co-formula (see (9j, Appendix A), the following formula for the free energy can 
then be obtained: up to an additive constant, 


A(z) = -p~ 1 ln(Z E J, 


(1.5) 


where Zs, = e 8^t x \_ z {dx), the submanifold £ z is defined by 


= {x = (xi ,..., x n ) G V | £(z) = z}, 

and 6^ X \- Z (dx) represents a measure with support £ z , such that 5^ x )_ z (dx)dz = dx 
(for further details on delta measures, we refer to [10] , Section 3.2.1). We assume 
henceforth that £ and V are such that Zs_ < oo, for all 2 £ M m . 

The idea of free energy biasing methods, such as the adaptive biasing force 
method (see 013 ) or the Wang Landau algorithm (see [14]), is that, if £ is well 
chosen, the dynamics associated to V — A o £ is less metastable than the dynamics 


f e^ ,3(V and |-Ad | denotes the Lebesgue 


associated to V. Indeed, from the definition of the free energy (1.4), for any compact 
subspace M C R m , the image of Z~ l e^~^ < 8 v ~ A °^ ( ' x ^\^ x ^ eM by £ is the uniform 

law 7 - 7 - 77 1 where Z = 

\M\ 

measure on A4. The uniform law is typically easier to sample than the original 
measure £ * fi. If the function £ is well chosen (i.e. if the dynamics in the direction 
orthogonal to £ is not too metastable), the free energy can be used as a biasing 
potential to accelerate the sampling of the dynamics (see [2j). The difficulty is 
of course that the free energy A is unknown and difficult to approximate using 


the original dynamics (1.2) because of metastability. Actually, in many practical 


cases, it is the quantity of interest that one would like to approximate by molecular 
dynamics simulations (see 0E3)- The principle of adaptive biasing methods is thus 
to approximate A (or its gradient) on the fly in order to bias the dynamics and to 


reduce the metastable features of the original dynamics (1.2). 


1.3. Adaptive biasing force method (ABF) 


In order to introduce the ABF method, we need a formula for the derivatives of 
A. The so called mean force X7A(z ), can be obtained from (1.5) as (see [TD], Sec¬ 
tion 3.2.2): 


V A(z)=f f{x)dpY. z , (1-6) 

J 

where dp^ z is the probability measure p conditioned to a fixed value 2 of the 
reaction coordinate: 


dgx z = Z^e-e v ^5 i{x) _ z (dx), 


(1.7) 
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and / is the so-called local mean force defined by 


fi = E G& Vfc.W - .1 ; d:v £ GfjVHj 

j-i \i=i 


( 1 . 8 ) 


where G = (Gij)ij = i i ..., rn , has components Gq ? - = • V£ 7 -. This can be rewritten 

in terms of conditional expectation as: for a random variable X with law ^ (defined 

by @), 

VA(z)=E(f(XMx) = z). (1.9) 

We are now in position to introduce the standard adaptive biasing force (ABF) 


technique, applied to the overdamped Langevin dynamics (1.2): 


dX t = -(\/V-J 2 F t° + V(W o £) ( X t )dt + ^2~F^dW tl 


i= 1 


( 1 . 10 ) 


F l(z) = E [fi( x t)\£( x t) = z], i = 1,..., m, 


where / is defined in (1.8). Compared with the original dynamics (1.2), two modi¬ 


fications have been made to obtain the ABF dynamics (1.10): 


(1) First and more importantly, the force FI o £V£i has been added to the 


i =1 


original force — XV. At time t, F t approximates VA defined in (1.6). 


(2) Second, a potential W of has been added. This is actually needed in the 
case when £ lives in an unbounded domain. In this case, a so-called con¬ 
fining potential W is introduced so that the law of f(X t ) admits a long¬ 
time limit Zf v 1 e~^ w ^dz (see Remark [ 2 ] at the end of Section 2.2), where 

Zw = / e~P w is assumed to be finite. When £ is living in a compact 

J Ran(£) 

subspace of K m , there is no need to introduce such a potential and the law of 
£(X t ) converges exponentially fast to the uniform law on the compact subspace 


(as explained in Section 1.2 and Section 2.2). Typically, W is zero in a chosen 
compact subspace M of R m and is harmonic outside M. For example, in di¬ 
mension two, suppose that £ = (£1,6) and M = [£ m m>£max] x If rn'i/n 5 (max] > 
then W can be defined as: 

2 2 

W (-Zl, Z 2 ) — ^ (Zi (max) T ^ (^x £mxxx) • (1-11) 


x= 1 


x= 1 


It is proven in [5] that, under appropriate assumptions, F t converges exponentially 


fast to VA. In addition, for well chosen the convergence to equilibrium for (1.10) 


is much quicker than for (1.2). This can be quantified using entropy estimates and 


Logarithmic Sobolev Inequalities, see [2j. 

Notice that even though F t converges to a gradient (VA), there is no reason why 
F t would be a gradient at time t. In this paper, we propose an alternative method, 
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where we approximate VA, at any time t, by a gradient denoted VA t . The gradient 
VA t is defined as the Helmholtz projection of F t . One could expect improvements 
compared to the original ABF method since the variance of VA t is then smaller 
than the variance of F t (since A t is a scalar function). Reducing the variance is 
important since the conditional expectation in ( flTOt is approximated by empirical 
averages in practice. 


1.4. Projected adaptive biasing force method (PABF) 

A natural algorithm to reconstruct A t from F t , consists in solving the following 
Poisson problem: 

A A t = divFt onA4, (1-12) 


with appropriate boundary conditions depending on the choice of £ and A4. More 
precisely, if £ is periodic and A4 is the torus T m , then we are working with periodic 
boundary conditions. If £ is with values in R m and A4 is a bounded subset of 
R m , then Neumann boundary conditions are needed (see Remark [8] at the end 
of Section 2.3.2). To solve this Poisson problem, standard methods such as finite 


difference methods, finite element methods, spectral methods or Fourier transforms 


can be used. Note that (1.121 is the Euler equation associated with the minimization 
problem: 


A t = argmin / \Vg-F t \ 
gem(M)/R JM 


(1.13) 


where H 1 (M)/M. = < g G H 1 (M) \ g = 0 > denotes the subspace of H 1 (M) of 


Ad- 


zero average functions. In view of (1.13), A t can be interpreted as the function 


such that its gradient is the closest to F t . Solving (1.12) amounts to computing the 


so-called Helmholtz-Hodge decomposition of the vector field F t as (see j6j, Section 3): 

F) = VAt + i?t, onAJ, (1-14) 


where R t is a divergence free vector field. 

Finally, the projected ABF dynamics we propose to study is the following non 
linear stochastic differential equation: 


( dX t = -V{V-A t oZ + Wo£)(X t )dt+VW :i dWt, 

< AA t = divFt onAJ, with appropriate boundary conditions, (1.15) 
l Ft( z ) = K [fi( x t)\C(X t ) =z\,i= 1 


Compared with the standard ABF dynamics (1.10), the only modification is that 
the mean force F t is replaced by VA t , which is meant to be an approximation of 
VA at time t. 

The main theoretical result of this paper is that A t converges exponentially 
fast to the free energy A in M (at least in a specific setting and for a slightly 
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modified version of (1.15), see Section [2] for more details). Moreover, we illustrate 
numerically this result on a typical example. From a numerical point of view, the 
interest of the method is that the variance of the projected estimated mean force 
(i.e. VA t ) is smaller than the variance of the estimated mean force (i.e. F t ). We 
observe numerically that this variance reduction enables a faster convergence to 
equilibrium for PABF compared with the original ABF. 

The paper is organized as follows. In Section [2j the longtime convergence of the 
projected ABF method is proven. Section [3] is devoted to a numerical illustration of 
the interest of the projected ABF compared to the standard ABF approach. Finally, 
the proofs of the results presented in Section [2] are provided in Section [4] 


2. Longtime convergence of the projected ABF method 

For the sake of simplicity, we assume in this Section that V = T n and that £ (x) = 
( x\,x 2 )■ Then £ lives in the compact space M. = T 2 and we therefore take W = 0. 
The free energy can be written as: 

A( Xl ,x 2 ) = -/3 _1 ln (Z S(a , ia , 2) ), (2.1) 

where = f e~ l3V< ' x ' ) dx 3 ...dx n and £( Xl>X2 ) = {xi_,x 2 } x T" -2 . The 

t '52(x 1 ,x 2 ) 

mean force becomes: 

Vd(ii,i 2 ) = f f(x)dnv (a:iiX2) , (2.2) 

where / = (fi,f 2 ) = (diF, d 2 V) and the conditional probability measure dpY. {xl } 
is: 

d Ms (xiiX2) = Z^ x2 e- pv dx 3 ...dx n . 

Finally, the vector field F t (x\,x 2 ) writes / /d/as ( , (t,.), or equivalently: 

F}(x!,x 2 ) =E(d z V(X t MX t ) = (x u x 2 j), *=1,2. 


2.1. Helmholtz projection 

In section |2.1.1| weighted Helmholtz-Hodge decomposition of F t is presented. In 
section [2.1.2| the associated minimization problem and projection operator are in¬ 
troduced. 

Let us first fix some notations. For x £ 1" and 1 < i < j < n, x? denotes 
the vector (xi,Xi+±, ... ,Xj ) and dx? denotes dxidxi+i ...dxj. Moreover, V x 2 , div x 2 
and A x 2 represent respectively the gradient, the divergence and the laplacian in 
dimension two for the first two variables (xi,x 2 ). Likewise, V x ™ = (d 3 , ...,d n ) T 
represents the gradient vector starting from the third variable of R". 
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2.1.1. Helmholtz decomposition 

The space T 2 is a bounded and connected space. For any smooth positive probability 
density function p : T 2 —> R, let us define the weighted Hilbert space: L 2 (T 2 ) = 
{/ : T 2 —► R, f T2 \f\ 2 <p < oo}. Let us also introduce the Hilbert space H v ( div; T 2 ) = 
{g £ L 2 (T 2 ) x L 2 (T 2 ), div x 2 (g) £ L 2 (T 2 )}. It is well-known that any vector field 
F t : T 2 —> R 2 £ iLi(div,T 2 ) can be written (see [Bj, Section 3 for example) as 
(Helmholtz decomposition): F t = A t + R tl where R t is a divergence free vector 
field. We will need a generalization of the standard Helmholtz decomposition to the 
weighted Hilbert spaces L 2 (T 2 ) and H v ( div;T 2 )): 

Ft<P = V x 2(H t )^j + R t , (2.3) 

s.t. div a , 2 (i? t ) = 0. This weighted Helmholtz decomposition is required to simplify 
calculations when studying the longtime convergence (see Remark [To] in Section [4. 1| 
for more details). Recall the space: iJ 1 (T 2 )/R = {g £ H 1 (T 2 ) | f j2 g = 0}. The 
function A t is then the solution to the following problem: 

[ V xl A t .V xl gy= [ FfV^gtp, Vg£H\ T 2 )/R, (2.4) 

J t 2 Jr 2 

which is the weak formulation of the Poisson problem: 


div x 2 (V x 2 A t ip(t ,.)) =div x 2 (F t ip(t,.)), (2.5) 

with periodic boundary conditions. Using standard arguments (Lax-Milgram the¬ 
orem), it is straight forward to check that (2.41 admits a unique solution A t £ 
H 1 ( T 2 )/R. 


2.1.2. Minimization problem and projection on a gradient 

Proposition 1 . Suppose that F t £ Ft v (div; T 2 ). Then for any smooth positive prob¬ 
ability density function ip, the equation ( |2.4[ ) is the Euler Lagrange equation asso¬ 
ciated with the following minimization problem: 


A t = min 

h£H 1 ( T 2 )/R. 


T 2 


\V x *h-F t \ 2 ip = 


mm 

heH 1 (T 2 )/R 


|V~? h — F t 


li 2 (T 2 )- 


( 2 . 6 ) 


Furthermore, A t belongs to H 2 ( T 2 ). 


Proof. Let us introduce the application I : H 1 ( T 2 )/R —> R + , defined by 

1(g) = 11 V x 2 h — F t 111,2 (T 2 )- It is easy to prove that I is a-convex and coercive, i.e. 

lim 1(g) = + 00 . Thus / admits a unique global minimum A t £ ft 1 (T 2 )/R. 
IMIhi-h-oo 
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Furthermore, \/e > 0, \/g £ Ft 1 ) T 2 ), 

I(A t +eg)= f (V x 2 (A t + eg) - F t \ 2 p 

J T 2 

= [ \^xf A t ~ F t \ 2 cp-2e [ (V x 2 A t — F t ) ■ V x 2 gip + e 2 j \V x 2 g\ 2 p 
J T 2 J T 2 J T 2 

= I(A t )-2ef (V x 2A t -F t )-V^gtp + e 2 j \X x 2g\ 2 p. 

J T 2 Jt 2 

(2.7) 

Since A t is the minimum of J, then I)A t + eg) > I)A t ), Ve > 0,Vg £ FI 1 ) T 2 ). 
By considering the asymptotic regime e —> 0 in the last equation, one thus obtains 
the equation (2.4): 

/ V x 2 A t • X x 2 gp = f F t ■ X x 2 gp, Mg £ J ff 1 (T 2 )/R. 

Jt 2 Jt 2 

This is the weak formulation of the following problem (2.5) in i?' 1 (T 2 )/R: 

div x 2 (V x 2 A t <p(t,.)) = div x 2 )F t p)t ,.)). 


Since ip is a smooth positive function, then 3d > 0, s.t. tp > S. Furthermore, since 
div x 2 )F t p)t,.)) £ L 2 ) T 2 ), thus A x 2 A t £ L 2 ) T 2 ). Therefore, using standard elliptic 
regularity results, A t £ Ft 2 ) T 2 ). □ 


For any positive probability density function p 1 the estimation vector field V x 2 A ( 
is the projection of F t onto a gradient. In the following, we will use the notation: 

V v (F t ) = V x? A t , (2.8) 

where the projection operator V v is a linear projection defined from H lp )div; T 2 ) to 
Ft 1 ) T 2 ) x H 1 ) T 2 ). Notice in particular that V v = V v . 


2.2. The projected ABF (PABF) method 

We will study the longtime convergence of the following PABF dynamics: 

f dX t = -V(V-A t o 0(X t )dt + ^JW'dWt, 

| V x 2 A t =V^i(F t ), (2.9) 

[ Ft(x u x 2 ) = E[diV(X t )\Z(X t ) = ( Xl ,x 2 )],i= 1,2, 

where is the linear projection defined by ( |2.8| ) and W t is a standard nN- 
dimensional Brownian motion. Thanks to the diffusion term \J2f3~ 1 dWt, X t admits 
a smooth density if) with respect to the the Lebesgue measure on T ra and ipt then 
denotes the marginal distribution of ip along £: 

ip*)t,x i,x 2 )= / ip)t,x)dx 3. 

-' s o 1 ,* 2 ) 


( 2 . 10 ) 
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The dynamics (2.9) is the PABF dynamics ( 1.15[ > with £(x) = (xi,x 2 ), W = 0 
and a weighted Helmholtz projection. The weight ip^ is introduced to simplify the 
convergence proof (see Remark [To] in Section El). 


Remark 1. If the law of X t is ip(t, x)dx then the law of £(X t ) is ip^(t, x\,x 2 )dxidx 2 


and the conditional distribution of X t given £(X t ) = (x\,x 2 ) is (see (1.7) for a sim¬ 


ilar formula when ip = Z. 


-l 


y v. 

(*1 .*2) 


-PV\. 


d[At,X 1,X2 — 


ip(t, x)dx 3 


Xi,X 2 ) ' 

Indeed, for any smooth functions / and g , 

®(f(Z(Xt))g(X t )) = f f(£(x))g(x)ip(t,x)dx 


( 2 . 11 ) 


0 


= / / / o ^gipdx^dxxdx-z 

J T 2 Js (xi!X2 ) 

r / v gipdx 3 

= / f{x\, x 2 ) --— V(x 1 ,x 2 )dx 1 dx 2 . 

J T 2 1p^{Xi,X 2 ) 


Let us now introduce the non linear partial differential equation (the so-called 
Fokker-Planck equation) which rules the evolution of the density ip(t, x) of X t solu¬ 
tion of (2.9): 


d t ip = div YjdiA t ofV&J ip + /3 l Xipj , for(t,a;) G [0,oo[xT", 

Vf > 0, div(VA t ?/^(f,.)) = div(F t ^(t,.)), in T 2 with periodic boundary conditions, 

[ diVipdx 3 

■,* = 1,2. 

( 2 . 12 ) 


Vf > 0, \/(x\,x 2 ) G T 2 , Fp(x\,x 2 ) = 


’-‘(x 1 ,x 2 ) 


ipt{x 1 ,x 2 ) 


The first equation of (2.12) rewrites: 

d t ip = di v[VVip + 1 V-0] - di({diA t )ip) - d 2 ((d 2 A t )ip). 


(2.13) 


Suppose that (ip,F t ,A t ) is a solution of ( 2 . 12 ) and let us introduce the expected 
long-time limits of ip, ip^ (defined by ( 2 . 10 l) and gt,xi,x 2 (defined by ( 2 . 11 )) respec¬ 
tively: 


( 1 ) j )ao=e -P{v-MQ. 

( 2 ) ip^ = 1 (uniform law); 

(3) 


HoO,Xi ,X2 


= z. 


-1 


-f> v da 


3 • 


(2.14) 
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Notice that the probability measure ip% 0 {x\,X 2 )dx\dx 2 is the image of the proba¬ 
bility measure ip 00 (x)dx by £ and that p.oo,xi,x 2 = Ms 


(*1»®2) 


defined in (1.7). Fur¬ 


thermore, we have that 


/ V^oo — 1? / 1 and V(^i, £ 2 ) € T , / dg 00jX ix 2 — 1* 

•'P Jt2 • / = ( . ll . 2) 


Remark 2. In the case when W ^ 0, the first equation of the Fokker-Plank prob¬ 


lem (2.12) becomes: 

d t ip = div (V (V — A t o £ — W o £) ip + /3 _1 V^) . 
The expected long-time limits of ip, ip^ and gt, Xl ,x 2 are respectively: 

(1) V’oo = Zrfe-K v - A °t- w °V-, 

( 2 ) 4 = 2 ^; 

(3) lloo,x 1: x 2 = Zy 1 e-^dx'i 


^(®1 .* 2 ) 

-pw 


3 I 


where Zw = I e 

J T 2 

2.3. Precise statements of the longtime convergence results 


In section 2.3.1 some well-known results on entropy techniques are presented. For a 
general introduction to logarithmic Sobolev inequalities, their properties and their 
relation to long-time behaviours of solutions to partial differential equations, we 


refer to ixicgna. Section |2. 3. 2| presents the main theorem of convergence. 


2.3.1. Entropy and Fisher information 

Define the relative entropy H(. |.) as follows: for any probability measures g and v 
such that g is absolutely continuous with respect to v (denoted g <C v), 

H{g\v) = J In dg- 

Abusing the notation, we will denote H(ig\ip) for H(tp(x)dx\ip(x)dx) in case of 
probability measures with densities. Let us recall the Csiszar-Kullback inequality 
(see [I]): 

\\g-v\\TV<y/2H(g\v), (2.15) 

where \\g— v\\tv = sup < / fd(g — v) > is the total variation norm of the signed 

ll/IU°°<i U J 

measure g — v. When both g and u have densities with respect to the Lebesque 
measure, \\g— v\\tv is simply the L 1 norm of the difference between the two densities. 
The entropy H(g\v) can be understood as a measure of how close g and v are. 
Now, let us define the Fisher information of g with respect to v. 


HhW) = J 


Ht) 


dg. 


(2.16) 
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The Wasserstein distance is another way to compare two probability measures 
H and v defined on a space E, 

W(ld,v) = i inf / ds(x, y) 2 dir(x, y), 
y frenOw Jsxs 

where the geodesic distance ds on E is defined as: Vr,|/ G S, 

dr (a;, 2/) = inf I J \w[t)\ 2 dt \ w G C^QO, 1], E), w( 0) = x, iu(l) = y 1, 


and v ) denotes the set of coupling probability measures, namely prob¬ 

ability measures on E x E such that their marginals are /i and v. W G 
Yl(^^)J^ x ^^(x)dTr{x,y) = and J Sx ^ip(y)dTT{x,y) = J^ifdv. 

Definition 1. We say that a probability measure v satisfies a logarithmic Sobolev 
inequality with constant p > 0 (denoted LSI(p)) if for all probability measure p 
such that p <C v , 

H(p\v) < M{p\v). 

2 P 

Definition 2. We say that a probability measure v satisfies a Talagrand inequality 
with constant p > 0 (denoted T(p)) if for all probability measure p such that p <C i', 

W(p,v) < ^-H(p\v). 

Remark 3. We implicitly assume in the latter definition, that the probability mea¬ 
sures have finite moments of order 2. This is the case for the probability measures 
used in this paper. 

The following lemma is proved in HU. Theorem 1: 


Lemma 1. If v satisfies LSI{p), then v satisfies Tip). 


Recall that X t solution to (2.9) has a density ip(t ,.). In the following, we denote 
the Total Entropy by 


E{t) = .)|V>oo) = / ln(V>Mx>)^, 

J T n 

the Macroscopic Entropy by 

E M {t) = .)|V4) = f In 

J T 2 

and the Microscopic Entropy by 

E m (t) = / e m {t,x 1 ,x2)i> i (t,xi,x 2 )dxidx2 , 

J M 


(2.17) 


(2.18) 


(2.19) 
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where e m (f, x±, X 2 ) = H(p t , Xl , X2 \poo,xi,x 2 )- The following result is straightforward 
to check: 

Lemma 2 . It holds, Vt > 0, 

E[t) = Em (t) + E m (t). 

Note that the Fisher information of pt,xi_,x 2 with respect to Poo,xi,x 2 can be 


written as (see (2.16)): 

I{St,Xl ,X 2 1 1^00,X\,X 2 ) — / | ^ X~ ^(^(t, ) | ~d[-lt,Xl r 


2.3.2. Convergence of the PABF dynamics (2.12) 


The following proposition shows that the density function satisfies a simple 
diffusion equation. 


Proposition 2 . Suppose that (il),F t ,A t ) is a smooth solution of (2.12). Then 
satisfies: 


( =/3 1 A x 2 ipS, in [0,oo[xT 2 , 

= V>o> on T 2 . 


( 2 . 20 ) 


Remark 4. If ^ = 0 at some points or is not smooth, then F at time 0 may 
not be well defined or I(ip^( 0, ■)/V’| 0 ) ma y be infinite. Since, by Proposition [ 2 ] ip£ 
satisfies a simple diffusion equation these difficulties disappear as soon as t > 0. 
Therefore, up to considering the problem for t > to > 0, we can suppose that i/jq is 

a smooth positive function. We also have that for all t > 0, ip^(t, .) > 0, / = 1 

J T 2 

and .) e C 00 ^ 2 ). 

Remark 5. In the case where W ^ 0, the probability density function satisfies 
the modified diffusion equation: 

d t V> ? = w). 

Here are two simple corollaries of Proposition [2] 

Corollary 1. There exists i 0 > 0 and J 0 > 0 (depending on such that 

Vf > to, -)|V4) < Io e ~ lj 1& ^ t - 


Corollary 2. The macroscopic entropy Em it), defined by (2.18), converges expo¬ 
nentially fast to zero: 

Vt > t 0 , E M (t ) < 

where Iq is the constant introduced in Corollary [7} 
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The assumptions we need to prove the longtime convergence of the biasing force 
VA t to the mean force VA are the following: 

[HI] V £ C 2 (T n ) and satisfies: 

3y > 0, V3 < j < n,Vx £ T n , max(|<9i9jF(a;)|, \dzdjV(x)\) < 7 . 


[H2] 14 is such that 3p 0 , N/(:z 7 ,:r 2 ) £ T , p , oo,xi,x 2 — Mh(xi,x 2 ) dehned by 
satisfies LSI(p). 

The main theorem is: 

Theorem 1. Let us assume [HI] and [H2[. The following properties then hold: 

(1) The microscopic entropy E m converges exponentially fast to zero: 

3C > 0,3A > 0,Vt > 0, y/E m (t) < Ce~ xt . (2.21) 

Furthermore, if p ^ An 2 , then A = /3 -1 min(p, 47 t 2 ) and 


C — \J E m ( 0 ) + -p— 


/3 1 |p — 47 t 2 | y 2 p 


—. If p = 4 - 7 T 2 , then for all A < /3 1 p, there 


exists a positive constant C such that (2.211 is satisfied. 


(2) y/E(t) and \\ijj(t ,.) - ^oo||i,i(T”) both converge exponentially fast to zero with 
rate A. 

(3) The biasing force X7 x 2A t converges to the mean force V x zA in the following 
sense: 

Vt> 0 , / \\/ x 2A t -V x 2A\ 2 ip i (t,x 1 ,x 2 )dxidx2<—E rn (t). ( 2 . 22 ) 

The proofs of the results presented in this section are provided in Section [4j 

Remark 6 . We would like to emphasize that our arguments hold under the as¬ 
sumption of existence of regular solutions. In particular, we suppose that the den¬ 
sity ip{t, .) is sufficiently regular so that the algebric manipulations in the proofs 
(see Section [4]) are valid. 

Remark 7. This remark is devoted to show how the rate of convergence of the 


dynamics (1.2 1 is improved thanks to PABF method. First of all, we mention a 


classical computation to get a rate of convergence for (1.2). Precisely, if one denotes 
p{t,.) the probability density function of X t satisfying (1.2), and ip 00 = Z~ 1 e~^ 


its longtime limit, then by standard computations (see for example 0) , one obtains: 

^-)l<Ax>) = —/3 _1 T(p(i, .)|y5oo). 

Therefore, if satisfies LSI(R), then one obtains the estimate 

3i? > 0, Vt > 0, flMi,.)l¥>oo) < H{p 0 \Too)e~ 2rlRt . (2.23) 


From (2.15), we obtain that ||<p(£,.) — < Poo\\L 1 n n ) converges exponentially fast to 
zero with rate /3~ 1 R. The constant R is known to be small if the metastable states 
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are separated by large energy barriers or if high probability regions for p are sepa¬ 
rated by large regions with small probability (namely p is a multimodal measure). 
Second, by Theorem[l] one can show that VA t converges exponentially fast to S/A 
in F 2 ('!/’| 0 (*i, x 2 )dx\dx 2 )—norm at rate A = /3 _1 min(p,47r 2 ). Indeed, since ip^ = 1, 


/ \S7A t -S7A\ 2 dx 1 dx 2 = [ \S7 A t — S/ A\ 2 ^^ XllX2 \ dxidx 2 

Jt 2 J T 2 ips{t,x l,X 2 ) 


< 


8y 2 


(l-e)p 




where e > 0 such that ip^(t,xi,x 2 ) > 1 — e (for more details refer to the proof of 
Corollary [l] in Section [4]). This result must be compared with ( |2.23[ ). More precisely, 
A is related to the rate of convergence 47 t 2 at the macroscopic level, for equation 
(2.201 satisfied by tp^, and the rate of convergence p at the microscopic level, com¬ 
ing from the logarithmic Sobolev inequalities satisfied by the conditional measures 
Moo,®i,a: 2 - Of course, p depends on the choice of the reaction coordinate. In our 
framework, we could state that a ’’good reaction coordinate” is such that p is as 
large as possible. Typically, for good choices of £, A 3> R, and the PABF dynam¬ 
ics converges to equilibrium much faster than the original dynamics ( |1.2| ). This is 
typically the case if the conditional measures p 0 o,x 1 ,x 2 are less multimodal than the 
original measure p. 


Remark 8. (Extension to other geometric setting) 

The results of Theorem [l] are easily generalized to the following setting: 

If V = £(x) = (xi,x 2 ) and M is a compact subspace of R n , then choose a 

confining potential W (defined in (1.11 1 ) such that Zw = f e~^ w < + 00 , Z/^e~^ w 
satisfies LSI{r*) (for some r* > 0) and W is convex potential, then Corollary [I] is 
satisfied with rate 2/? _1 (r* — e), for any e € (0, r*) (refer to Corollary 1 in [9] for 
further details). In this case, Neumann boundary conditions are needed to solve the 


Poisson problem (2.5): 


div(VA t i/> c (f,.)) = d\v{F t ^{t ,.)) 
dA t 


dn 


= Ft.n 


in M, 


on dM, 


(2.24) 


where n denotes the unit normal outward to A4. The convergence rate A of The¬ 
orem |T] becomes /3 -1 min(p, r* — e). Neumann boundary conditions come from the 
minimization problem ( |2.6[ ) associated to the Euler-Lagrange equation. The numer¬ 
ical applications in Section [3] are performed in this setting. 


Remark 9. (Extension to more general reaction coordinates) 

In this section, we have chosen £(xi,..., x n ) = (aq, x 2 ). The results can be extended 
to the following settings: 


(1) In dimension one, the Helmholtz projection has obviously no sense. If T) = T n 
and £(a;) = X\, then F t converges to A', which is a derivative of a periodic 
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function and thus / A' = 0. Since 

J T J T 

therefore take a new approximation A t = F t — 

J T 

at any time t. The convergence results of this section can be extended to this 
setting, to show that A' t converges exponentially fast to A !. 

(2) More generally, for a reaction coordinate with values in T m , the convergence 
results presented in this paper still hold under the following orthogonality con¬ 
dition: 


F t is not necessary equal to zero, one can 
F t , which approximates A! 


Vi ^ 3 , V& • Vi', = 0. 


(2.25) 


In the case when (2.251 does not hold, it is possible to resort to the following 


trick used for example in metadynamics (refer to [3J [8]). The idea is to introduce 
an additional variable z of dimension m, and an extended potential (x, z) = 
V(x) + ^\z — £(x)\ 2 , where k is a penalty constant. The reaction coordinate is 
then chosen as ^ me ta(x,z) = z, so that the associated free energy is: 


At(z) = -r 1 In/ e~P v ^ x ’ z) dx, 
Jv 


which converges to A(z) (defined in (1.51) when n goes to infinity. The extended 
PABF dynamics can be written as: 

dX t = - |vk(l t ) + k - Z it tMi(X t )^j dt + VW^dWt, 

dZ t = n(£(Xt) - VE t (Z t )) dt + VW 1 dW t , 

XE t = V^ meta ( G t ), 

[G t (z)=E(aX t )\Z t = z), 

where Wt is a m —dimensional Brownian motion independent of W*. The results 
of Theorem [l] apply to this extended PABF dynamics. 


3. Numerical experiments 

3.1. Presentation of the model 

We consider a system composed of N particles (qi)o<i<N-i in a two-dimensional 
periodic box of side length L. Among these particles, three particles (numbered 0, 1 
and 2 in the following) are designated to form a trimer, while the others are solvent 
particles. In this model, a trimer is a molecule composed of three identical particles 
linked together by two bonds (see Figure [I]). 

3.1.1. Potential functions 

All particles, except the three particles forming the trimer, interact through the 
purely repulsive WCA pair potential, which is the Lennard-Jones (LJ) potential 
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truncated at the LJ potential minimum: 


VwcA(d) 



iid < do, 

if d > do, 


where d denotes the distance between two particles, £ and <r are two positive pa¬ 
rameters and do = 2 1 / 6 cr. 

A particle of the solvent and a particle of the trimer also interact through the 
potential Vwca ■ The interaction potential between two particles of the trimer (< 70 /q\ 
or Q 1 /Q 2 ) is a double-well potential (see Figure [ 2 ]): 


V s (d) = h 



(d — d\ — w) 2 1 2 


(3.1) 


where d\, h and w are positive parameters. 

The potential Vs has two energy minima. The first one, at d = d\, corresponds to 
the compact bond. The second one, at d = d\ + 2 oj, corresponds to the stretched 
bond. The height of the energy barrier separating the two states is h. 

In addition, the interaction potential between qo and <72 is a Lennard-Jones 
potential: 


V LJ (d) = 4e' 




where s' and a' are two positive parameters. 

Finally, the three particles of the trimer also interact through the following 
potential function on the angle 9 formed by the vectors qiq 0 and <71 <72: 

Vo 0 (9) = y (cos( 6 ») - cos( 6» 0 )) 2 , 

where 0q is the equilibrium angle and kg is the angular stiffness. Figure |T| presents 
a schematic view of the system. 



Fig. 1. Trimer (go* <71? 92 )- Left: compact state; Center: mixed state; Right: stretched state. 
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The total energy of the system is therefore, for q G (LT) 2Ar : 

2 TV —1 

V(q) = £ VwcA(\qi - qj |) + EE VwcA{\qi - qj |) 

3<T<j<TV— 1 2—0 j— 3 

1 

+ E V s(\qi - *+l|) + - g2|) + Vfl o (0). 

i=0 



Fig. 2. Double-well potential (|3.1|), with d\ = 2 1 / 6 , tu = 2 and h = 2. 


3.1.2. Reaction coordinate and physical parameters 

The reaction coordinate describes the transition from compact to stretched state in 
each bond. It is the normalised bond length of each bond of the trimer molecule. 
More precisely, the reaction coordinate is £ = (£1,62)) with £1 ( q ) = ' qo ~^ j~ da an d 
£ 2 ( 9 ) = ^-^T da ■ For i = 1,2, the value ^ = 0 refers to the compact state (i.e. 
d = d 0 ) and the value = 1 corresponds to the stretched state (i.e. d = d 0 + 2 w). 

We apply ABF and PABF dynamics to the trimer problem described above. 
The inverse temperature is /3 = 1, we use N = 100 particles (N — 3 solvent particles 
and the trimer) and the box side length is L = 15. The parameters describing the 
WCA and the Lennard-Jones interactions are set to a = 1, e = 1 , o' = 1, e' = 0.1, 
do = 2 1 / 6 , c?i = 2 1 / 6 and the additional parameters for the trimer are ui = 2 and 
h = 2. The parameters describing the angle potential are: do such that 003 ( 00 ) = 1/3 
and kg = 1 (we refer to m, Section 10.4.2, for the choice of such parameters). The 
initial condition on the trimer is as follows: Both bonds qoqi and qiq 2 are in compact 
state, which means that the distance between qo and q\ and the distance between 
qi and q 2 are equal to do- Moreover, the initial bond angle is 0 q- 
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3.1.3. Numerical methods and numerical parameters 

Standard and projected ABF methods are used with N rep n cas = 100 replicas of the 
system evolving according to the overdamped Langevin dynamics discretized with 
a time-step At = 2.5 x 10 -4 . The reaction coordinate space of interest is taken 

of the foim A4 = [£mm;£max] ^ [£mim ^max\i where £,rnin ~ 0• - and t^rnax — 1-2. 

A4 is discretized into Ni, ins x N^. ins = 50 x 50 = 2500 bins of equal sizes and 
S = 5 X = S y = = 0.028 denotes the size of each bin along both axes. 

To implement the ABF and PABF method, one needs to approximate FI (x, y) = 
E[fi(x,y)\t;(x,y) = (£ 1 ( 2 ;, y), £ 2 ( 2 , ?/))], i = 1,2. The mean force F t is estimated in 
each bin as a combination of plain trajectorial averages and averages over replicas. 
It is calculated at each time as an average of the local mean force in the bin over 
the total number of visits in this bin. More precisely, at time t and for l = 1,2, the 
value of the mean force in the ( i,j) th bin is: 


Fibs) 


Nreplicas 

E E fl (Qk,t'')l{indx(€(q k t /))=(i,j)} 

t'<.t k—1 

N rep i icas 

E E ('))=(*.!)} 

t'<t k—l 


(3.2) 


where qk,t. denotes the position ( Xk,yk ) at time t, f = (/i,/ 2 ) is defined in ( 1 . 8 ) 
and indx(^(qk,t’)) denotes the number of the bin where £(qk,t:) lives, i.e. 


indx{^(q)) 



£l (q) £min 


£,min 


_i 

+ 

-1 

_1 


\/q G M. 


If the the components of the index function (i.e. indx) are either equal to —1 or to 
Nbins, it means that we are outside A4 and then the confining potential is non zero, 
and defined as: 


2 

W(£(qk,t)) = ^ — £max) + 1 { £ 4 (q k _ t ) < f mi „ } (£i ( Qk , t ) — in) ] ■ 

i—1 


To construct the PABF method, the solution to the following Poisson problem with 
Neumann boundary conditions are approximated: 


AA = divF 

¥ = F 

on 


n 


in M = [£„ 
on dM. 


n £>max\ ^ Antin’) ^ s max\’> 


(3.3) 


where n denotes the unit normal outward to M. We use for simplicity in the nu¬ 
merical experiments the standard Helmholtz problem, without the weight ijfi. Prob¬ 


lem (3.3) is solved using finite element method of type Q 1 on the quadrilateral mesh 
defined above, with nodes ( Xi,yj ), where Xi = £, m in + id and yj = £ m in + j8, for 
i,j = 0, Nbins- The space M is thus discretized into Nt = N^ ins squares, with 
N s = (Nbins + l ) 2 nodes. The associated variational formulation is the following: 
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f Find A e H l {M)/M. such that 

| [ VA-Vu= / F ■ Vv, Vue H 1 (M)/R. 

[ J M JM 

3.2. Comparison of the methods 

In this section, we compare results obtained with three different simulations: without 
ABF, with ABF and with projected ABF (PABF). First, it is observed numerically 
that both ABF methods overcome metastable states. Second, it is illustrated how 
PABF method reduces the variance of the estimated mean force compared to ABF 
method. As a consequence of this variance reduction, we observe that the conver¬ 
gence of VA t to VA with the PABF method is faster than the convergence of F t 
to VA with the ABF method. 


3.2.1. Metastability 

To illustrate numerically the fact that ABF methods improve the sampling for 
metastable processes, we observe the variation, as a function of time, of the two 
metastable distances (i.e. the distance between go and ■ and the distance between 
gi and 52 )- On Figures [3] and |dj the distance between q\ and <72 is plotted as a 
function of time for three dynamics: without ABF, ABF and PABF. 

Both ABF methods allow to switch faster between the compact and stretched 
bond and thus to better explore the set of configurations. Without adding the 
biasing term, the system remains trapped in a neighborhood of the first potential 
minimum (i.e. do ~ 1.12) region for 20 units of time at least (see Figure [3]), while 
when the biasing term is added in the dynamics, many jumps between the two local 
minima are observed (see Figure [4]). 



time 


Fig. 3. Without ABF. 
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Fig. 4. Left: ABF; Right: PABF. 


3.2.2. Variance reduction 

Since we use Monte-Carlo methods to approximate VA t , the variance is an impor¬ 
tant quantity to assess the quality of the result. The following general proposition 
shows that projection reduces the variance. 


Proposition 3. Let F be a random function from T 2 into R 2 and belongs to 
i?(div,T 2 ), and define V = V\ (i.e. without weight) the projection on gradient 
vector fields defined in Section \ 2.1.2 Then, the variance ofV(F) is smaller than 
the variance of F: 


Var (V{F)) < / Var(F), 


IM 


IM 


where, for any vector field F , Var(.F') = E(|F| 2 )—E(|F |) 2 and\F\ being the Euclidian 
norm. 


Proof. Let F be a random vector field of H{ div,T 2 ). Let us introduce V(F) £ 
H l { T 2 ) x H 1 ( T 2 ) its projection. Notice that by the linearity of the projection 
V(E(F)) = E (V{F)). By definition of V{F), one gets: 


[ {F-V(F))-X7h = 0,Vh£H 1 (T 2 ). 

Jr 2 

Therefore, using Pythagoras and the fact that V{F) is agradient, 

/ \F\ 2 = [ \F-V(F)\ 2 + [ \V(F)\ 2 

Jr 2 J t 2 J t 2 


and 


[ |F-E(F)| 2 =/ \F-E(F)-V{F-E{F))\ 2 + I \V(F-E(F))\ 2 . 

J T 2 JT 2 JT 2 

Using the linearity of V, we thus obtain 

/ Var(F) = / Var(P - V{F)) + / Var (V(F)), 

J t 2 Jr 2 jt 2 

which concludes the proof. 


□ 
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We illustrate the improvement of the projected method in terms of the variances 
of the biasing forces, by comparing Var(VA t ) = Var(9i A t ) + Var((9 2 A t ) (for the 
PABF method) and Var(F t ) = Var(F t 1 ) + Var(.F 2 ) (for the ABF method). Figure [5] 
shows that the variance for the projected ABF method is smaller than for the 
standard ABF method. 

We have Nuns x Nb ins = 2500 variable for each term (i.e. d\A t , 82 At, F} and 
i 7 ) 2 ). The variances are computed using 20 independent realizations as follows: 


1 50 1 20 50 / 20 

Var Vt) = E ™ £ F }' k &> yrf ~ 2500 £ ( 9n E F t 


2500 ^ 20 

1 , 3=1 k —1 


ij =1 


20 




*:=1 


Note that four averages are involved in this formula: an average with respect to the 
space variable, an average over the 20 Monte-Carlo realizations, an average over 
replicas and a trajectorial average (the last two averages are more explicit in (3.2)). 
Notice that since the variance of the biasing force is smaller with PABF, one may 


expect better convergence in time results. This will be investigated in Section 3.2.3 
and Section 13.2.41 



Fig. 5. Variances as a function of time. 


3.2.3. Free energy error 

We now present, the variation, as a function of time, of the normalized L 2 — distance 
between the real free energy and the estimated one, in both cases: ABF and PABF 
methods. As can be seen in Figure [Gj in both methods, the error decreases as time 
increases. Moreover, this error is always smaller for the projected ABF method than 
for the ABF method. 
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Fig. 6. Free energy error as a function of time. 


3.2.4. Distribution 

Another way to illustrate that the projected ABF method converges faster than the 
standard ABF method is to plot the density function ip* as a function of time (see 
Figure [t]|TT|). It is illustrated that, as time increases, the probability of visiting all 
bins (of the reaction coordinate space A4) increases. 

It is observed that, for the projected ABF method, the state where both bonds 
are stretched is visited earlier (at time 5) than for the standard ABF method (at 
time 20). The convergence to uniform law along (^1,^2) is faster with the projected 
ABF method. 




Fig. 7. At time 0.025. Left: f il>^{x\,X 2 )dx 2 ', Right: f if>€(xi,X 2 )dxi. 
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Fig. 8. At time 5. Left: f / i/j*(xi,X 2 )dx 2 ', Right: f if*(xi,X 2 )dxi. 




Fig. 9. At time 10. Left: f X 2 )dx 2 ‘, Right: f (x±, X 2 )dx±. 



Fig. 10. At time 20. Left: f if€(xi, X 2 )dx 2 ', Right: f if€(x±,X 2 )dx\. 


4. Proofs 

The proofs are inspired from [5]. One may assume that /? = 1 up to the following 
change of variable: t = /3 _1 £, ?/>(£, x) = x), V(x) = /3V(x), Recall, that we work 

in V = T n , M = T 2 and Vx = (xi, x n ) G T n , £(x) = (xi,X2). 
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Fig. 11. At time 25. Left: f X2)dx2] Right: f X 2 )dx\. 


4.1. Proof of Proposition ^ ] 

Let g : T 2 —>■ K. be a function in H 1 ( T 2 ). 


d f df 

— / ip i gdxidx 2 = — / i^g ° Cda;™ 
dt Jj2 dt Jjn 

= [ div[(w- + 

Jt- i=1 

/. 2 2 

= - / Et( w - E ^ ° ev&w + v#v^p o ^ 

jTa j =1 i=l 

= -T I 

i=l “'T" 

+ E l diA t ° £,i’d z g ° (,dx 

Jt" 


Applying Fubini, it holds: 

7 /. 2 


dt L ^ gdxidx 2 = -^2 


/ T 2 


= 1 V ' T2 v ' S (^ 1,®2) 
2 


[c^Vd/' + 9i'0]dx3 9ip(a;i, x 2 )dx±dx 2 


+E 


— 1 "'^' 2 
2 


diA t {xi,X2)^dx^dig{xi,X2)dxidx2 

^ 2 
= -E / Ft'*P i d i g(x 1 ,x 2 )dx 1 dx 2 - V' / di^d i g(x 1 ,x 2 )dx 1 dx 2 

i=1 i=l 7t2 

+ E / diAtfxi^^ip^digixxjX^dxidx^ 

i=i ^ 

= / A^g(x 1 ,x 2 )dx 1 dx 2 , 

JT 2 
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where we used (2.4) with ip = ^(t ,.). This is the weak formulation of: 

d t V = A V, on [0, oo[xT 2 . 

Remark 10. The reason why we consider the weighted Helmholtz decomposi¬ 
tion (2.3) with ip = in the PABF dynamics (2.9) instead of the standard 

one (1.14) is precisely to obtain this simple diffusion equation on the function i[>^. 
This is will also be useful in the proof of Lemma [6] below. 


4.2. Proof of Corollary [7] 

Let p = ipt and p oo = = 1. It is known that Vi > 0 and V(xi,X 2 ) € T 2 , p 

satisfies: 


dtp = A^. 

Moreover (See Remark [4]), it is assumed that and is such that 


(4.1) 



landed,.) > 0. 


Let us show that Vi > 0,V/c > 0, || p(t ,.) — l|| J Ff*(x 2 ) < ||0(O,.) — l||//fc(T 2 ) e Sn2t . 
First, we prove that p converges to 1 in L 2 ( T 2 ), 


1 d_ f 

2 djt J r ^2 


\cP~l \ 2 



dtP(P - 1) 
AM-1) 

I V0V(0-1) 

T 2 



<-4t r 2 / \<j> -1| 2 , 
J T 2 


where we have used the Poincare-Wirtinger inequality on the torus T 2 , applied to 
p: for any function f £ H l { T 2 ), 



We therefore obtain, || p(t ,.) - 1|| 2 2(T 2) < ||</>(0, ■) - l|| 2 2 (T 2 )e 

Second, we prove that dip converges to 0 in L 2 (T 2 ). For i = 1,2, dip satis- 
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fies (4.11: dt(di<j)) = A x 2 (di<j>), with periodic boundary conditions. As above, 

f l^l 2 = f dtidifydi# 

z at j T 2 j T 2 

= I A (di<j))di(j) 

J T 2 

= -[ |V(5^)| 2 . 

J T 2 


[ \di<j>\ 2 < -4?r 2 [ f di(f> — [ di(j> 

J T 2 JT 2 V JT 2 


Using again Poincare-Wirtinger inequality on 

1 d 

2 dt _ _ 

= -4tt 2 [ \di(f>\ 2 . 

J T 2 

Where we used f j2 di<f> = 0, since (j> is periodic on T 2 . Therefore, it holds 

l|di<HMIl!2 ( T 2 ) < ll^^(0,-)lli 2 (T2) e“ 87r2 ‘- 

Third, one can prove by induction that all higher derivatives of <j> converge 
exponentially fast to 0, with rate 87 t 2 and the following estimation is then proven: 

Vt > 0,Vfc > 0, •) ~ l||^ffc( T 2 ) < ||<^(0,.) — l||_f/'i(T 2 ) e 

As H k { T 2 ) ^ L°°(T 2 ), Vfc > 1, then 3c > 0, 

110 - lllioo < c||0 - l||^ fe < ce" 8 - 2 *. 

Therefore, Ve > 0. 3 1 0 > 0, Va; £ T 2 , Vf > t 0 , 6(t, x) > 1 — £. Finally, Vi > t 0 

f f ^ 2 <l|V, 2 0(O,)||| 2(T2 )_ 87r2t 

/T 2 0 1 — £ Jf 2 Xl 




1 — £ 


4 . 3 . Proof of Corollary ^2 ] 

We have that V4o = 1 satisfies LSI(r), for some r > 0 (see Chapter 3, Section 3 in 
PP). Referring to Proposition [ 2 ] and since is a probability density function, one 
gets: 

j t E M = j d t (V> £ InC0 £ )) 

= [ d t ijfi In(^ £ ) + f d t ^ 

= f ln(f/>^) 

J T 2 

= -/ |V x 2ln(^)| 2 ^ 

Jt 2 

= —2rE M - 
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Therefore, Em converges exponentially fast to zero. Referring to Corollary [I] and 
since Em converges to zero, we have that for any t > to, 

/ OO 7 f- oo 

—E M (s)ds = -Jt I(V\*l>io)ds 


— 8 tts 


-~ Io J t e 

— p -8ti - 2 t 

8tt 2 


ds 


which yields the desired estimation. 


4.4. Proof of Theorem [T] 

To prove our main result, several intermediate lemmas are needed. 

Lemma 3. Vi > 0, V(xi,x 2 ) £ T 2 and for i = 1,2, we have: 

(F[-d z A)(x i,x 2 ) = ^J 5 i ln(V’/V’oo)^a;3j {x l ,x 2 )-{d 1 \n{V/Voo)) (* 1 .^ 2 ) 


Proof. 

L 


di ln('il>/'ip 00 )^dx% - di ln(^VV4 


di hi^f-dxo - [ di In (ijj 00 )-j-dx% - d z ln(V> £ ) + di ln(^) 


1 

V JY 


(®1 > x 2 ) 


diipdxo + f (diV - V&diA o £)^-dx% - di ln(^ s ) 


^ + F ; _ a ( A - *** 






= Ft - diA. 


□ 


Lemma 4. Suppose that [HI] and [H2] hold, then for all t > 0, for all (x’i, x 2 ) £ T 2 
and for i = 1,2, we have: 


|i r t *(xi,x 2 ) - 5iA(xi,x 2 )| < yy-e m (i,xi,x 2 ). 

Proof. For any coupling measure 7r £ J^(/Ut, Xl , X2 ,Moo,xi,a: 2 ) defined on £( xliX2 ) x 
S( xllX2 ), it holds: 


1 ^ - = 


(diV(x) — diV{x'))-K{dx,dx') 


= ||V x? ^P|| L oo / / 


^(*1 ,a! 2 ) X ^(a:i .* 2 ) 


ds (xiiI2) (%, x') 2 n(dx, dx'). 
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Taking now the infimum over all 7 r G H(M*, .|(xi,x 2 )),/i(oo, .|(*i) ^ 2 ))) and using 
Lemma [I] we obtain 

I FI - diA\ < 7 W(p(t, .\(x 1 ,x 2 )),^(oo, .\{x 1 ,x 2 ))) 


< l\j ~H(^(t, .\(xi,x 2 )), pt(oo, .{(xx^^)) 


= 7 \/ (xi,x 2 )). 


Lemma 5. Suppose that [H2] holds, then for all t > 0, 

E m (t)<^-[ |V a ,jln(V'(t,.)/V’oo)| 2 '0- 
2 P J T" 

Proof. Using [H2], 

E m = e m ^dxidx 2 
J T 2 

= / .|(a:i,a;2))|Ai(oo, ,|(a;i,a;2)))-!/' £ rfa:irfa;2 

Jt 2 


< 


f 1 

/T 2 , 


|V E? ln(V>(t,.)/V>oo)| 2 rfa:3 ^ 


’ lpt(t,Xi,X 2 ) 


dxidx 2 


= 7T I |V x? ln(V»(t, O/V'oo)! 2 ^"- 

2 P J T" 


□ 


□ 


Lemma 6. It holds for all t > 0, 

[ (diA t - Fl)[di ln(?/;/V>oo)]t/> + f (d 2 A t - F?)[d 2 ln(V>Mx>)]V> < 0. 

JT n JT 71 

Proof. Using Fubini, 

f (diAt - Fl)[di ln(‘if>/'i/j 00 )\ip + [ (d 2 A t - F?)[d 2 ln(ip/ip^ip 
JT n J T n 

= f (diAt-F}) f [5iln(^/^ 00 )]V , + [ ( d 2 A t -F.f) [ [0 2 hi(^Mx>)]'*/’- 

Jt2 • / =Ci.. a > 

For the first term, we have 

(di In if)ip- j (d\ In il>oo)ip 
d ®(®1 ,a:2) 

= + [ di(V — A)ip 

= (d\ In ip£)ijj£ + F/ip^ — d\ AipZ. 


[0i \n(ip/ip^ip = 


x 2 ) 
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Similarly, we have 



[d 2 ln(V>/V’oo)]V’ 


(d 2 ln?/^)V^ + — d 2 Axj}^. 


Therefore, one gets 


/ (diA t - Fl)[di ln('0/V>oo)]^ + / (c^t - F 2 )[d 2 

Jjn Jjri 

= [ {diA t - Ft){di \nV)V + [ {d 2 A t -F?)(d 2 \nV)V 

J T 2 Jt 2 

- [ {hAt - F/)V - [ (d 2 A t - i?)V 

Jt 2 Jt 2 


which concludes the assertion since the first line is equal to zero (by (2.4) with 
ip = ijj^{t,.)) and the second line is non positive. Again the weighted Helmholtz 
decomposition helps in simplifying terms. □ 


Proof of Theorem [l} 


Now we will prove the exponentially convergence of E m {t). Recall (2.13): 


dti> = div(W+ Vi/0 - di{{d\A t )il>) - d 2 {{d 2 A t )ip), 


which is equivalent to 


d t ip = div ( V’ooV ( ) ) + 9i[(9i A - diA t )i/j\ + d 2 [(d 2 A - d 2 A t ) 


Using (2.17), (2.18) and (4.1), one obtains 


dE 

dt 


dE M 

dt 


|Vln(V>/V’oo)| 2 t/’ + / {diAt-d 1 A)[d 1 \n{il)/tl> 00 )]il> 

n Jjn 

(d 2 A t - d 2 A)[d 2 ln^/V’oo )]ip, 


— f |V ln(il>£)\ 2 V■ 
J T 2 
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Using then Lemma [2] and Lemma [3j one gets 
dE m dE dEM 

dt dt dt 

= ~ f iVln^/V’oo)! 2 ^ + f (diA t -d 1 A)d 1 \n(ip/-ip 00 )ip 

Jj n Jjn 

+ f {d 2 A t - d 2 A)d 2 ln(^i/V’oo)V’ + / |<9i lm/> 5 | 2 ^ + f \d 2 ln^ 5 | 2 -i/> £ 

J T" J T 2 J T 2 

= -[ iVln^/V’oo)! 2 ^ 

,/T 71 

+ [ (diA t - Fl)[di lii('0/V>oo)]^ + f {Ft - d 1 A)[di\n(il)/ip 00 )]'il) 

Jjn Jjn 

+ f {d 2 A t - F 2 )[<9 2 ln^/VtxOlV’ + [ ( F t ~ d 2 A)[d 2 ln(WV’oo)]V> 

+ f |9i + f \d 2 

J T 2 Jt 2 

Lemma [6] then yields 
dE, 


dt 


<- [ |Vln(^oo)|V 

JT 71 

+ [ {Ft - d x A)di ln^/^oo)^ + [ ( F f - d 2 A)d 2 ]n{ip/ijj 00 )ip 

Jjn Jjn 


/ T 2 


|9i ln^| 2 ^ + / |c? 2 ln^| 2 ^. 


/T 2 


Using lemma [3] and Fubini, one then obtains 

iVln^/^oo)! 2 ^ 


/T 2 


di ln(WV’oo) 

L, (xi,aT 2 ) 


/T 2 


±_ 

1p£ 

f 9 2 \n(t/j/ipoo)^ 

^( x l ! a3 2 ) ^ 


5i \n{'ij)/ip 00 )ip - d i ln(^)<9i ln(V>Mx>)^ 


-‘(xi ,x 2 ) 


a 2 ln(V’/V’oo)^ - d 2 ln{^)d 2 ln{ip 
J ^‘(x 1 ,x 2 ) J T 71 


/T 2 


|c?! ln?/^| 2 ^ + / |<9 2 ln'i/^| 2 'i/^ 


/T 2 


<-/ (Vhi^/V’cxOlV 


/t 2 


,x 2 ) 


<9i ln(^/V>oo)^ 


/T 2 


d 2 ln(ip/'if) 00 )ip 


Ll ( x l .^2) 


-jj ~ f d\ \n(ip^)di ln(V’/'0oo)V ) 

^ JT 71 

[ d 2 ln(^)d 2 \n{ip/'ip 00 )ip 


1 


/T 2 


|9i ln^| 2 ^ + / |3 2 ln^| 2 V^- 


/T 2 
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Applying Cauchy-Schwarz on the first terms of the second and third lines, we obtain 

dEm < - f |V x n ln(^oo)| 2 ^ 

JT n 

- [ d\ ln(?/^)9i ln^/V’oo)^ — f d 2 ln(tp^)d 2 ln(ip/i/j 00 )ip 

Jjn Jjn 

+ f \d\ lnV\ 2 V + [ \d 2 lnV\ 2 V 

J T 2 J T 2 


dt 


< - 


|V*» ln(V’/V’oo)| 2 V’ - [ di ln(V> £ ) [ <9i ln^/V’oo)-^ - <9i ln^ £ 




- f d 2 ln(^> 5 ) 
JT 2 


[ d 2 ln(V’/^oo) 77 - d 2 In V 


ip. 


Applying Lemma [5j Lemma [3J Cauchy-Schwarz, Lemma [4] and Corollary [T] 

dEr, 


dt 


<-2 pE m + J j \di ln(^)| 2 ^«, / [ -e m (t,(xi,x 2 ))^ 
y Jt 2 y Jt 2 P 


/T 2 


|S 2 ln(^)| 2 V’ £ A 


/t 2 P 


&m {t, (xi,x 2 ))i>t 


< -2pE m + 2j ] J^E rn ^jJ^ |V a ,2 ln(-^ £ )| 2 -!/' £ 

< -2P-E m + 2 ^y - p E m\lliV/^to) 

< —2pE m + 2jyJ~^~ErnVhie 4,r 


Finally we obtain 


— jE m (t) < -p^/E m {t) + 74 / f-e 


0 47T 2 t 


dt 


2 P 


First, if p 7 ^ 4 - 7 T 2 , using Gronwall inequality, one obtains 
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\J E m (t) < \J E m ( 0) e 



which leads the desired estimation ( |2.21[ ). 

Using this above convergence, Corollary [2] and Lemma [2] it is then easy to see 
that E converges exponentially fast to zero. Using ( 2. 15| ) , one obtains the conver¬ 
gence of to ipoo since: 

\\4> - ^ooIIlht") < = V2E. 

The second point of the theorem is checked. Finally, we are now in position to prove 
the last point of Theorem [I] Using (2.6) and Lemma |4j 

||VA t - VA||^ (t2) < 2||VA t — -ft||+ 2|| F t - VA||^2^ T 2) 

<4||F t -VA||^ (T2) 

7 2 

< 8— E m . 
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